
/**** Step 1:  Create moments, normal***/

clear all
set more off


use "$path\Intermediary Data\localshock_dataforReg_sp95_v2.dta",clear
drop  j j3 local_shock10 m_dprix
drop gasoline_rotterdam_euro dprix1 


drop if prix_ht==.
*** Step 1 ***

gen t=1
bysort id: egen nobs_stations=sum(t)

keep if nobs_stations>500


gen minus_number_obs=-nobs_stations
sort minus_number_obs id


**** Step 2 : Iterate ****

egen id_select=group(minus_number_obs id)
drop minus_number_obs
drop if id_select==.

bysort id (daten): gen dp=log(prix_ht[_n])-log(prix_ht[_n-1])
bysort id (daten): gen df=log(prix_ht[_n])-log(prix_ht[_n-1]) != 0
gen dp_nonnul=dp if dp!=0
gen dum_dp=(dp!=0)
replace dum_dp=. if dp==.

gen dp_nonnulp=dp if dp>0
gen dp_nonnuln=dp if dp<0

drop dp
sum dp_nonnul*, d


keep id daten dp_nonnul dp_nonnulp dp_nonnuln dum_dp id_select

preserve
gen t=1
collapse  (mean) freq=dum_dp  (sum) nobs=t, by(id)
save "$path\Intermediary Data\moments_EG_sp95", replace
restore

drop if dp_nonnul==.
drop dum_dp
sum id_select
local max=r(max)
forvalue j = 1/`max'  {

preserve
keep if id_select==`j'
capture drop num_obs
gen num_obs=_n


winsor2 dp_nonnul , replace  cuts(1 99) trim

qui replace dp_nonnulp=. if dp_nonnul==.
qui replace dp_nonnuln=. if dp_nonnul==.

gen abs_dp=abs(dp_nonnul)
gen t=(dp_nonnul!=.)
gen tp=(dp_nonnulp!=.)
gen tn=(dp_nonnuln!=.)

egen ndp=sum(t)
egen ndpp=sum(tp)
egen ndpn=sum(tn)


egen kurtosis=kurt(dp_nonnul)
egen mean=mean(dp_nonnul)
egen skewness=skew(dp_nonnul)
egen sd=sd(dp_nonnul)
egen meanp=mean(dp_nonnulp)
egen meann=mean(dp_nonnuln)
egen mean_abs=mean(abs_dp)

egen p25=pctile(dp_nonnul), p(25)
egen p50=pctile(dp_nonnul), p(50)
egen p75=pctile(dp_nonnul), p(75)

egen p25a=pctile(abs_dp), p(25)
egen p50a=pctile(abs_dp), p(50)
egen p75a=pctile(abs_dp), p(75)

egen p25p=pctile(dp_nonnulp), p(25)
egen p50p=pctile(dp_nonnulp), p(50)
egen p75p=pctile(dp_nonnulp), p(75)

egen p25n=pctile(dp_nonnuln), p(25)
egen p50n=pctile(dp_nonnuln), p(50)
egen p75n=pctile(dp_nonnuln), p(75)



bysort id (daten): gen lag1dp_nonnul=dp_nonnul[_n-1]
bysort id (daten): gen lag2dp_nonnul=dp_nonnul[_n-2]
bysort id (daten): gen lag3dp_nonnul=dp_nonnul[_n-3]
bysort id (daten): gen lag4dp_nonnul=dp_nonnul[_n-4]
bysort id (daten): gen lag5dp_nonnul=dp_nonnul[_n-5]

qui reg dp_nonnul lag1dp_nonnul lag2dp_nonnul lag3dp_nonnul lag4dp_nonnul lag5dp_nonnul
gen autocor=_b[lag1dp_nonnul]
****

qui duplicates drop id_select, force
keep id kurtosis skewness  mean* p* sd ndp ndpp ndpn auto

if `j'==1 {
save "$path\Intermediary Data\moments_alt_sp95.dta", replace
}
if `j'>1 {
append using "$path\Intermediary Data\moments_alt_sp95.dta"
qui save "$path\Intermediary Data\moments_alt_sp95.dta", replace
}

restore
drop if id_select==`j'

display "End of iteration number" `j'
}


**** Step 2: corrected kurtosis ****

use "$path\Intermediary Data\localshock_dataforReg_sp95_v2.dta",clear
drop  j j3 local_shock10 m_dprix
drop gasoline_rotterdam_euro dprix1 


drop if prix_ht==.
*** Step 1 ***

gen t=1
bysort id: egen nobs_stations=sum(t)

keep if nobs_stations>500



gen minus_number_obs=-nobs_stations
sort minus_number_obs id


**** Step 2 : Iterate ****

egen id_select=group(minus_number_obs id)
drop minus_number_obs
drop if id_select==.

bysort id (daten): gen dp=log(prix_ht[_n])-log(prix_ht[_n-1])
bysort id (daten): gen df=log(prix_ht[_n])-log(prix_ht[_n-1]) != 0
gen dp_nonnul=dp if dp!=0
gen dum_dp=(dp!=0)
replace dum_dp=. if dp==.

gen dp_nonnulp=dp if dp>0
gen dp_nonnuln=dp if dp<0

drop dp
sum dp_nonnul*, d




keep id daten dp_nonnul dp_nonnulp dp_nonnuln dum_dp id_select



drop if dp_nonnul==.
drop if dp_nonnul==0

drop dum_dp dp_nonnulp dp_nonnuln
sum id_select
local max=r(max)
forvalue j = 1/`max'  {

preserve
keep if id_select==`j'
capture drop num_obs



winsor2 dp_nonnul , replace  cuts(1 99) trim

qui drop if dp_nonnul==.
sort daten
gen n=1
egen t=group(daten)

xtset n t
keep n t dp_nonnul id  id_select 


*1. Demean dlog 
egen mdlogp=mean(dp_nonnul)
gen dprix_m=dp_nonnul-mdlogp

**2. dlogp2 dlogp4
gen dlogp4=dprix_m^4

gen dlogp2=dprix_m^2
gen mis=(dlogp2!=.)

sort n t
foreach i of numlist 1/25{
xtset n t
gen l`i'_dlogp2=l`i'.dlogp2
gen mis`i'=(l`i'_dlogp2!=.)
gen prod`i'=l`i'_dlogp2*dlogp2
egen s_p`i'=sum(prod`i')
egen s_mis`i'=sum(mis`i')
}

egen N_4=mean(dlogp4)
 
foreach i of numlist 5(5)25{

	gen D_sum`i'=0
	gen T_sum`i'=0

	foreach k of numlist 1/`i'{
qui	replace D_sum`i'=D_sum`i'+s_p`k' 
qui	replace T_sum`i'=T_sum`i'+s_mis`k'
	}

	gen denom`i'=D_sum`i'/T_sum`i'
	gen k_`i'=N_4/denom`i'
}




egen kurtosis=kurt(dp_nonnul)

qui duplicates drop id_select, force
keep id kurtosis k_*

if `j'==1 {
save "$path\Intermediary Data\kurto_hetero_sp95.dta", replace
}
if `j'>1 {
append using "$path\Intermediary Data\kurto_hetero_sp95.dta"
qui save "$path\Intermediary Data\kurto_hetero_sp95.dta", replace
}

restore
drop if id_select==`j'
display "End of iteration number" `j'
}
